Introduction
Millennial-length networks of annually resolved past-climate proxies derived from tree rings, ice cores, and historical records are pivotal in assessing the role of temperature and precipitation variability (LaMarche, 1974; Williams and Wigley, 1983; Graumlich, 1991; Bradley and Jones, 1992; Briffa et al., 1992; Graumlich, 1993; Scuderi, 1993; Stine, 1994; Hughes and Graumlich, 1996; Mann et al., 1998; Mann et al., 1999; Barber et al., 2000; Crowley, 2000; Stahle et al., 2000). Especially important contributors to the time-series data are tree-ring records from high-elevation, long-lived conifers in western North America, and high-latitude conifers (LaMarche et al., 1984; Briffa et al., 1992; Graumlich, 1993; Graybill and Idso, 1993; Mann et al., 1998; Mann et al., 1999; Briffa and Osborn, 1999; Barber et al., 2000). It is especially important to include high-elevation tree-ring series in large-scale reconstructions because 20th-century trends in high-elevation growth rates, including CO2 fertilization (Graybill and Idso, 1993; Jacoby and D'Arrigo, 1997), have been attributed as significant factors in addition to temperature. Some concerns stem from the widespread use of strip-bark trees in dendroclimatological reconstructions, which might show enhanced growth rates in the 20th century (Jacoby and D'Arrigo 1997). Strip-bark formation (i.e., trees with partial cambial dieback) is a morphological characteristic commonly displayed by trees growing at high elevations (Fig. 1). Strip-bark trees might show enhanced CO2 fertilization effects due to increased carbon allocation to a smaller associated surface area in the active cambium (Graybill and Idso, 1993). Increased ambient CO2 might also reduce stomatal conductance in strip-bark trees, resulting in higher water-use efficiency and subsequent fertilization effects (LaMarche et al., 1984; Graumlich, 1991). Analyses and interpretations of CO2 fertilization effects on strip-bark trees are controversial. Notably absent is acknowledgment of abiotic environmental correlates and spatial variability in strip-bark tree distribution.
A high-elevation whitebark pine (Pinus albicaulis Engelm.) site in the Greater Yellowstone Ecosystem was discovered to have very old strip-bark and entire-bark trees growing in proximity on similar slopes. This data set, containing 27 pairs of strip-bark and entire-bark trees, allowed for direct comparisons of the influence of tree growth form on tree-ring growth rates over time and space. These data can expand basic understanding of strip-bark tree autecology through temporally and spatially explicit analyses of where, when, and how tree locations and growth rates relate to the abiotic environment. While the field of landscape ecology has been founded in the spirit of answering such questions (Risser et al., 1984; Urban et al., 1987), those techniques are not in wide use by paleoecologists.
We conducted a spatial analysis of the distribution pattern of the strip-bark trees to determine the distribution as compared to a random expectation. After finding spatial structure to the strip-bark tree distribution, we assessed the strength of environmental correlates to presence or absence of strip-bark trees using a spatially explicit regression. We also examined the growth patterns of the strip-bark trees as compared to their entire-bark neighbors and paid particular attention to the role of an abiotic proxy for soil moisture. This analysis illuminates the complexity of climate-growth relationships within a species and highlights the importance of considering tree physiognomy and microsite variation in developing climatically sensitive tree-ring chronologies.
Materials and Methods
STUDY AREA
The study area is a 450,000-m2 (45 ha) rectangle in Montana's Tobacco Root Mountains, centered on 45°34′20″N, 112°20′25″W (Fig. 2). It is a steep and uniformly forested open-canopy subalpine site approximately 200 m below alpine treeline, with an average elevation of 2850 m. The site is a monospecific stand of whitebark pine. A census of the area was conducted for strip-bark trees. Each of the 27 strip-bark trees was located and mapped with a GPS receiver and differentially corrected to within 5 m horizontal accuracy. To analyze the growth rates of the strip-bark trees while accounting for abiotic factors, an entire-bark companion tree (defined as being <10 cm difference in diameter at breast height) was located for each strip-bark tree. The average distance between a strip-bark tree and its entire-bark companion was 4.6 m (range 0.97–6.7 m).
TIME-SERIES DATA
Increment cores were taken at breast height (1.37 m) from all the trees in the study. Two cores were taken from each strip-bark tree: one from the living cambium and one from the dead side. The core was taken from the dead side of the tree in an attempt to coarsely date the initiation of strip-bark morphology. When possible, cores were taken perpendicular to the slope to minimize the effects of slope pressure on wood formation. Although not quantified, the aspect of living cambium was not obviously systematic (i.e., strip-barking occurred at multiple aspects irrespective of slope facet). A core was taken from each companion tree as well. All of the cores were cross-dated and measured using standard methods (Fritts, 1976; Stokes and Smiley, 1996). Raw ring widths were used, as opposed to normalized widths, to highlight absolute changes in growth rates within and between trees. Ring widths are most appropriate in this and similar dendroecological studies (cf. Fritts and Swetnam, 1989) when ecological and/or physiological processes at the tree or stand level are being investigated. Analyzing ring widths allowed for the variance between strip-bark and entire-bark trees to be emphasized and not equally scaled over the analysis period, as a ring-width index would have done.
A yearly mean growth rate was calculated for the strip-bark and entire-bark trees from 1750 to 1999 A.D. The year 1750 was selected as a cutoff date so that we could compare the 125 years before and after the two types of trees began to diverge. Limiting the analysis to 1750 also reduced the age-related growth trends in the individual series, and therefore the need to detrend. To test for significant differences in growth before and after the industrial revolution, the chronologies were divided into two periods: the early period spanned 1750–1874 A.D., and the later period spanned 1875–1999 A.D. After both mean-ring-width series were found to have significant temporal autocorrelation at lag of 2 yr, autoregressive moving average models (ARMA) were fit to remove that autocorrelation (Box et al., 1994). A paired t-test was used on the residuals from the second-order ARMA model to test whether a difference existed between the mean of the strip-bark and entire-bark ring-width series (Zar, 1996).
To assess the difference between the strip-bark and entire-bark trees, the mean growth difference (MGRD) between the strip-bark and entire-bark tree rings was calculated as:
where SBi and EBi were ring widths in mm at year i for strip-bark and entire-bark pairs, and n was the number of years in the chronology. This metric gave a measure of the difference between strip-bark and companion entire-bark trees over period n (1750 to 1999 A.D.).
ENVIRONMENTAL DATA
A 10-m raster digital elevation model (DEM) was obtained from the U.S. Geological Survey ( http://www.gisdatadepot.com/dem). Slope, topographic convergence, and potential relative radiation were calculated from the DEM. A topographic convergence index (TCI) measured the tendency of water to collect on the landscape. It was calculated as:
where α is the upslope contributing area from the watershed and β is the local slope for each cell ij (Moore et al., 1991; Urban et al., 2000). The index takes on high values in coves or streambeds and low values on ridgetops and other drained areas.
A potential relative radiation (PRR) index measured the amount of sunlight that a particular raster element received given the shading that its neighborhood provided (cf. Allen and Walsh, 1996). It was calculated as:
where Ak was the sun angle and Bk was the sun azimuth measured at time step k, m was the number of time steps, βij was the local slope, and σij was the local aspect at cell ij. The total number of cells was represented by n. This index was calculated for the entire viewshed every 15 minutes over the growing season and then rescaled to 8 bits (0 to 255). (Ken Pierce and colleagues at the Landscape Ecology Lab at Duke University developed the PRR algorithm \[Pierce et al., unpublished manuscript\]. A macro was written in Arc Macro Language to perform this task in ArcGIS [ESRI 2003] and is available from Bunn.)
Each strip-bark tree was described by a data frame of abiotic variables (slope, TCI, and PRR). For pairs of strip-bark and entire-bark trees where it was possible to ascertain the date of strip-barking (n = 15) MGRD was added to the data frame. A series of 27 random points was generated in the study area and described by the abiotic variables.
POINT PATTERN
The distribution pattern of the strip-bark trees was described using Ripley's K function (Ripley, 1976, 1981). Ripley's K describes the cumulative frequency distribution of observations at a given point-to-point distance. Simulation envelopes for the pattern were plotted using 1000 Monte Carlo simulations of a Poisson process and provided a test against complete spatial randomness (CSR) (Manly, 1997). A buffer was invoked to restrict the analysis to distances of half the smallest dimension of the study area (Cressie, 1993). Transformed as L̂, the expected value under CSR is a straight line (equal to the distances). Values of L̂ greater than the distance indicate clumps, while values less than the distance indicate regularity (Ripley 1976).
MANTEL'S TEST
Initial inspection of the residuals of a logistic regression of the environmental variables on tree presence indicated significant autocorrelation (not shown). Results from logistic regressions with nonindependent residuals are not reliable, and we therefore chose an alternative analytic framework that is robust to spatial structure. Using a series of Mantel's tests, we identified variables that had the highest degree of correlation with the observed strip-bark tree location as compared to random points. The Mantel's test assessed the degree of correlation between strip-bark tree location and underlying environmental variables (slope, PRR, and TCI) while taking into account the relative spatial location of these variables and the intercorrelation among other variables. In a Mantel's test, the variables are dissimilarity matrices, and the test measures the degree of pairwise similarity between samples (Mantel, 1967; Legendre and Fortin, 1989; Leduc et al., 1992). This measurement examines whether samples that are similar in environmental variables also are similar in terms of strip-bark tree occurrence. A Mantel's test also considers geographic location as a predictor variable, indicating whether samples that are close together in space are similar in other variables.
A powerful form of Mantel's test is a partial regression using species dissimilarity (or in this case, presence or absence of strip-bark trees), abiotic dissimilarity, and geographic dissimilarity (distance) (Smouse et al., 1986; Legendre and Fortin, 1989; Fortin and Gurevitch, 1993). This form indicates how much variability in species composition and other variables are explained by abiotic factors. It also indicates spatially structured residual variability in species composition after removing the effects of the environmental variables.
Five distance matrices were constructed. The space, slope, PRR, and TCI matrices were constructed using Euclidean distance. The tree-distance matrix was constructed as a group contrast matrix, where the input variable from the raw data array was scored as 0 if strip-bark trees were absent from the sample location and 1 if they were present. That is, similar sites had a contrast of 0, and dissimilar sites had a contrast value of 1 (Schick and Urban, 2000). Because the elements of a distance matrix are not independent, a test of significance was evaluated via permutation (Legendre and Fortin, 1989; Manly, 1997). Mantel's tests were also performed on the same suite of environmental variables regressed against the strip-bark and entire-bark trees' relative growth rates (MGRD). Again, five Euclidean dissimilarity matrices were used.
Results
TIME SERIES
Although cores were taken from the dead side of all 27 strip-bark trees, we could only cross-date a subset of 15 because of rot on the dead sides of the remaining trees. The date of strip bark for these trees was extremely variable, ranging from 1621 to 1940 A.D. (Fig. 3). The date of the strip-barking event is an estimate because of erosion of the dead wood exposed to the elements (see Discussion).
The time-series plots of the 15 datable strip-bark chronologies and their neighbor entire-bark companions show nearly identical growth rates and coherent trends until about 1875 A.D. (Fig. 4). At that time, the mean strip-bark widths increase but largely follow the interannual to decadal pattern of the entire-bark chronologies. Box plots demonstrated obvious differences before and after 1874 A.D. (Fig. 5). There is no significant difference between the strip-bark and entire-bark chronologies from 1750 to 1874 A.D. as measured by a paired t-test on the residuals of the ARIMA model (t = −1.65, p-value = 0.11). The strip-bark trees, however, do show a significant increase over the entire-bark trees from 1875 to 1999 A.D. (t = 2.34, p-value = 0.02).
DISTRIBUTION PATTERN
Although tree locations within the study area were evenly distributed, the strip-bark tree locations were shown to be more clustered than random using Ripley's K. At distances up to 175 m, the values of L̂ were substantially greater than the simulation envelopes from 1000 Monte Carlo simulations of CSR (Fig. 6). This difference indicates a spatial pattern in strip-bark tree distribution that is not evident in the entire-bark trees.
MANTEL'S TEST
The regression of slope, PRR, and TCI on tree presence indicated that TCI was a significant predictor of tree presence (Table 1, Fig. 7). When controlling for spatial structure and the other variables, TCI was the only important patterning variable. A significant pure spatial residual was also present in this data, indicating an unspecified variable controlling the distribution of the trees or the random points.
The regression of slope, PRR, and TCI on MGRD demonstrated that TCI was a significant predictor of differences in tree growth (Table 2, Fig. 8). When controlling for spatial structure and the other variables, TCI (a proxy for soil moisture) was the most important patterning variable.
Discussion
These data provided a first look at how strip-bark vs. entire-bark tree distribution within a site can help disentangle the factors leading to strip-bark tree presence and growth. The causes of strip-bark formation and the subsequent interpretation of these trees' growth are complicated and lead to questions about tree physiology, climate, and atmospheric chemistry. We have added another layer of important factors to this question by showing that the spatial distribution and growth rates of these ancient trees were not random and were correlated to an abiotic proxy for soil moisture. The analysis of tree distribution data, when coupled with the time-series analysis, allowed us to explore mechanisms underlying the growth patterns. The spatial analyses lent support to the idea that abiotic processes and spatial heterogeneity are important.
Although increment cores were extracted from the dead side of all 27 of the strip-bark trees, it was possible to get an approximate date of the strip-bark event on 15 stems. This was due to external weathering on the dead wood exposed to the elements and internal rot making the precise dating of the strip-barking event using dendrochronological methods problematic. Furthermore, on these 15 stems the date of the strip-barking event is a point estimate given that cambial dieback likely occurs over decades, spreading around the bole. Given this, we limited the time-series analysis to these 15 stems to conservatively interpret the growth rates between strip-bark and entire-bark trees.
The strip-bark time series showed a consistent pattern of postindustrial growth rates similar to that found in the White (Graybill and Idso, 1993; Tang et al., 1999) and Sierra Nevada Mountains (Graybill, 1987; Bunn, unpublished data) of California. The tree-ring patterns in Figure 4 are the most extensive pairing of strip-bark and entire-bark trees published to date and supported the concept that strip-bark tree growth was unlike entire-bark tree growth after about 1875 A.D. The multiple dates of cambial dieback (Fig. 3) precluded an age-related trend in strip-bark tree growth pattern. The synchronicity of the departure in time and across space suggested a driving mechanism at regional or larger scales. There was no evidence, however, that a subcontinental climate forcing was responsible for a growth departure in 1875. Instead, the departure was concurrent with the onset of the industrial revolution and the associated rise in CO2. The time series did not show a steady linear increase in growth rates. Rather, it exhibited several rapid increases in growth over 2 to 5 yr, followed by variation for 1 or 2 decades around the new mean. This pattern was consistent with many ideas about complexity and threshold responses in ecosystems (Holling, 1973; Scheffer et al., 2001).
When we initially examined the tree locations in Figure 2, a strong spatial pattern was evident. That pattern suggested clustering at short distances and perhaps regularity at longer distances. The plot of L̂ also indicated this pattern by showing clustering at all scales less than 175 m (Fig. 6). The slight dip at 250 m indicated that the L̂ measurement might exhibit regularity at longer scales, were the data able to support this analysis. The effect of the edge of the study site, however, precluded reasonable analysis at distances greater than 250 m (Cressie, 1993).
The most important interpretation of the Ripley's K analysis was that the strip-bark tree distribution was not random. Absence of a complete stem map at the site precluded a quantitative analysis of the entire stand, but in-field examination demonstrated that the study area was continuously and homogenously forested, with trees occurring evenly throughout. The presence of entire-bark trees throughout the area and strip-bark trees within clusters indicated that something was driving the strip-bark tree locations. It was conceivable that a genetic predisposition to strip-bark growth form existed, and tree clustering was occurring as a function of seed dispersal. It was more likely, however, that an abiotic proxy was controlling strip-bark tree distribution.
The Mantel's test on presence (Table 1, Fig. 7) of strip-bark trees supported the hypothesis that strip-bark locations were controlled by abiotic factors. TCI, an abiotic proxy for soil moisture, was the significant environmental variable in this data set. It was logical to hypothesize that this sample of strip-bark trees was clustering on wet places because the trees were individuals that would have died in other physical niches but continued to grow in the most favorable locations remaining. Forest simulation models show tree distribution retreating into more favorable niche spaces in times of climatic stress (Brubaker, 1986).
If these strip-bark trees were clustering in relatively wet areas, then the important issues remaining were to understand strip-bark tree growth rates in terms of the abiotic environment and possible mechanisms for CO2 influence. It has been hypothesized that CO2 fertilization might increase with elevation as the partial pressure of CO2 decreases. There is small evidence that some alpine forbs are CO2 limited (Körner, 1999) but little evidence that increases in elevation directly affect tree growth (Graumlich, 1991; Jacoby and D'Arrigo, 1997). A more likely pathway for CO2 to affect tree growth is through modification of water-use efficiency via the negative response of stomatal conductance (Strain and Cure, 1985). Increases in water-use efficiency with increasing CO2 have been demonstrated experimentally (Idso et al., 1985; Mooney et al., 1991), but a field experiment designed to test for differences in water-use efficiency between strip-bark and entire-bark trees found none (Tang et al., 1999). Our data could not give a definitive answer to water-use efficiency question. They did, however, point to some interesting trends involving spatial heterogeneity.
Tree location in terms of species composition and relative growth rates has been shown to be controlled by abiotic variables (Oberhuber and Kofler, 2000). The Mantel's test on the MGRD differences (Table 2, Fig. 8) showed that the difference in growth rates increased with the proxy for soil moisture (increasing TCI). This conclusion was consistent with a much simpler hypothesis that xylem growth might be amplified geometrically as the proportion of crown and root biomass to active xylem increases. However, such a hypothesis does not explain the departure of the strip-bark growth rates after about 1875 A.D. that are seen in multiple records and of trees that strip-barked long before 1875 A.D. An intriguing possibility is that the water-use efficiency is coupled to CO2 increase and that the response of tree growth is strongly nonlinear, with a step change occurring as a function of CO2 around 1875 A.D. However, our data cannot explicitly explore this prospect.
CONCLUSIONS
These data demonstrated the utility of exploring the complex relationship between temporal and spatial data in the context of seeking to understand multiple interacting factors governing 20th-century growth anomalies. The wealth of climate information in dendrochronology is a critical data source. It is also apparent that the physical template is an important patterning agent in controlling time-series data.
Strip-bark growth rates departed from their entire-bark companions around 1875 A.D. Furthermore, the strip-bark tree distribution was not random but was correlated with an abiotic proxy related to soil moisture. Finally, the differences in growth rates between strip-bark trees and their entire-bark companions increased with likely increased soil moisture, as represented by TCI.
Although these data showed the importance of accounting for microsite variability when considering tree growth, they were not appropriate to rule definitively on the appropriateness of strip-bark trees as recorders of climate. For instance, calibration of high-elevation chronologies to lower-elevation chronologies has proven effective at dealing with tree-ring data (Mann et al., 1999), and relatively simple standardization techniques (e.g., splines) are also effective. Rather than speak to the robustness of tree-ring chronologies, we hope that these data will spark the debate on the relationship between physiognomy, spatial heterogeneity, and growth rates.
Acknowledgments
Jeremy Littell and Sean Hill at the Big Sky Institute's Tree-Ring Lab assisted in collecting and analyzing the tree-ring data. Dean Urban at the Landscape Ecology Lab at Duke University and Dr. Sara Goslee of the USDA-ARS Jornada Experimental Range at New Mexico State University provided the ECODIST S-Plus library used to execute the Mantel's tests. It is available from Dr. Urban ( http://www.env.duke.edu/landscape). Rob Schick and Cathy Zabinski provided helpful insights during the preparation of this manuscript. The Canon National Parks Science Scholars Program (Bunn), the Montana Space Grant Consortium (Bunn and Graumlich), and the Undergraduate Scholar's Program at Montana State University (Bellante) provided funding for this work.
References Cited
FIGURE 1.
One of the 27 strip-bark whitebark pines used in the study. The death of the cambium on the right side is complete, while the left side of the tree continues to grow
FIGURE 2.
Map of the study area in the Tobacco Root Mountains of Montana. Strip-bark tree locations are shown on a 20-m contour map
FIGURE 3.
The approximate date of the strip-barking event is shown for 15 of the 27 trees used in the study. Although cores were taken from the dead sides of all trees, only the 15 shown were possible to date using dendrochronological methods
FIGURE 4.
Time-series plot of the mean growth increment (in mm) for the 15 datable strip-bark trees and their entire-bark companions
FIGURE 5.
Box and whisker plots of the mean yearly growth of the 15 strip-bark trees and their entire-bark trees companions from 1750 to 1874 and from 1875 to 1999 A.D. The strip-bark and entire-bark chronologies were significantly different for the latter period only
FIGURE 6.
Plot of L̂ values marked as triangles at a given distance in meters. The lines are simulation envelopes that are the result of 1000 Monte Carlo simulations
FIGURE 7.
Results from Mantel's test for presence/absence data presented as a path diagram. Arrows indicate a significant effect. Magnitude of arrow represents strength of Mantel coefficient. The arrows to the left of the ecological variables indicate simple correlations, while the arrows to the right indicate partial correlations (compare to Table 1)
FIGURE 8.
Results from Mantel's test for growth data presented as a path diagram. Arrows indicate a significant effect. Magnitude of arrow represents strength of Mantel coefficient. The arrows to the left of the ecological variables indicate simple correlations, while the arrows to the right indicate partial correlations (compare to Table 2)
Table 1
Results of simple and partial Mantel's tests for presence of strip-barked trees.a TCI was the most important variable once spatial autocorrelation was taken into account (n = 54)
Table 2
Results of simple and partial Mantel's tests for relative growth rates of strip-barked trees.a Again, TCI was the most important variable once autocorrelation was taken into account (n = 15)
